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SECTION 1 


INTRODUCTION 

The dynamics of a helicopter blade in forward flight are usually 
described by a system of differential equations with periodic coefficients. 

A growing acceptance of hingeless helicopter blades for conventional 
helicopters flying at relatively high forward flight speeds has intensified 
the need for fundamental research on the aeroelastic stability of such 
systems. 

Studies dealing with the effect of forward flight (or periodic 
coefficients) have been primarily devoted to the study of flapping insta- 
bility at high advance ratios (Refs. 1 through 8). A limited number of 
studies dealing with the effect of periodic coefficients on coupled flap- 
lag (Refs. 9 and 10) or coupled flap-lag-pitch (Ref. 11) motion were also 
conducted. The case of coupled flap-lag motion has been, somewhat incon- 
clusively, investigated by Hall (Ref. 10) using multivariable Floquet 
theory. The same problem was also considered by Friedmann and Tong (Ref. 9) 
but the treatment was limited to low advance ratios (y < 0.3). Coupled, 
linearized, flap-lag-torsion motion has been investigated by Crlmi 
(Ref. 11) using a modified Hill method. In both cases (Refs. 10 and 11) 
only a limited number of numerical results were obtained and the physical 
mechanism of the aeroelastic instabilities has not been clearly identified. 

In particular, the degree of freedom which is responsible for the instability 
was not identified and the results for forward flight were not compared 
with those for hover. 
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Recent investigation of the aeroelastlc stability of hingeless blades 
in hover (Ref. 12) indicated that the aeroelastlc stability boundaries are 
quite sensitive to the number of degrees of freedom employed in the analy- 
sis. Therefore it is important to determine how the flapping behavior of 
a blade at high advance ratios is modified by the inclusion of the lag 
degree of freedom. This important problem, which has not received adequate 
treatment before, is one of the main topics of the present study. 

The mathematical methods used in previous studies dealing with the 
effects of forward flight were: 

(a) The rectangular ripple method (Ref, 1) 

(b) Analog computer simulation (Refs, 3 and 4) 

(c) Various variations of Hill's method (Refs, 2 and 11) 

(d) Multivariable Floquet-Llapunov theory (Refs. 6, 7 and 10) 

(e) Perturbation method in multiple time scales (Refs, 8 and 9) 

The mathematical method employed in the present study is the Floquet- 
Liapunov theorem, and the transition matrix is evaluated by two separate 
methods : 

1, Direct niimerlcal integration using a fourth order Runge- 
Kutta method 

2. A new and computationally efficient method developed by Hsu 
(Refs. 23 through 25) which is a multivariable extension of the 
rectangular ripple method. 

It is also shown that careful use of these methods enables one to circum- 
vent problems associated with identifying the results encountered in 
previous studies (Ref. 10). 

In addition, a new and convenient approximation for the reversed flow 
region is developed. This approximation is believed to be adequate for 
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most blade stability analyses. Finally, the effects of various Important 
parameters such as collective pitch setting, structural damping, droop and 
precone on the stability associated with forward flight are investigated. 



SECTION 2 

THE EQUATIONS OF MOTION 

2.1 Basic Assijfflptlons 

The geometry of the problem is shown in Fig. 1. The following basic 
assumptions were used in deriving the equations of motion: (a) The blade 

is cantilevered at the rotor hub. It can have an angle of droop Bjj at the 
root. In addition, the feathering axis can be preconed by an angle 0p. 

The angles and Bp are small, (b) The blade is torsionally rigid and can 
bend In two directions normal to the elastic axis, (c) The deflections of 
the blade are moderately small so that terms of 0(£jj ) can be neglected com- 
pared to one. (d) Moderately large deflections have only a small effect on 
the tension due to elastic effects since one of its ends is free, therefore 
a linear treatment of the elastic restoring forces is adequate, (e) A two- 
dimensional aerodynamic strip theory is used with C(k) * 1 and apparent mass 
effects are neglected, (f) Stall and compressibility effects are neglected, 
(g) Reversed flow is Included using an approximate model for reversed flow 
described in Appendix D. (h) A single blade analysis is performed, where the 
blade is attached to an aircraft with infinite mass. 

In deriving the equations of motion, an x, y, x coordinate system 
(Fig, 1) rotating with the shaft of the helicopter and attached to the 
blade is used. 

2.2 Derivation of the Equations of Motion 

The equations of dynamic equilibrium of a blade undergoing only 
bending in flap and lag can be taken from Friedmann (Ref. 13): 
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where and P^ Include the aerodynamic and inertia loads distributed 
along the span of the rotor blade. 

Denoting the elastic coupling effect by and defined in 

Appendix E, Eqs, (E,3) and (E.4), and assuming that e * 0, Eqs. (2.1) 

A 


can be rewritten as 

.2 
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( 2 . 2 ) 


(2.3) 


o o 

The loading terms in the x, y and z directions with nonlinearities up 
to and including second order terms in displacements can be written as 
(Reference 13) 


3T „2 

= -m 52 
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(x^ + e^^ + j/J - 2vJ 


P * L -m 52' 

y y 


V - (v) + 2uj - gg^ 52 


(2.4) 

(2.5) 
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(2.6) 


2 A 

-m w - ggp « 

where terms marked by are small^ and have been neglected. The last 
terms in Eqs* (2.5 and 2.6) represent viscous type structural damping. 
The boundary conditions for a cantilevered helicopter rotor blade 


are 


3w (o,t) 9^w (l,t) (l,t) 

“e ' -3F TT- n ® 


3x 




(2.7) 


av^(o,t) 3^v^(i,t) (i,t) 

^ ^- TT ■ = — 0 


9x 


3x*“ 


o o 

The total displacements in the x,y,z directions are given by 
(Ref. 13)^^ 

2 v2l 


u = -w^ (e„ + ** T 

e p d 2 


“.•V -5 


dx^ (2.8) 


~ (2.9) 

w = + x^ (3p + Bjj) (2.10) 

where the last term in Eq. (2.8) represents the shortening effect due to 
bending under the assumption of an inextensible blade. 

Applying Galerkin's method on the equations of motion, the spatial 
variable is eliminated. This is achieved by representing the elastic part 
of the displacements v and w by a sum of normal modes in flap and lag 
respectively. 


with assumption C p.4 
^^with sin 0 “ 6 
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( 2 . 12 ) 


”e ■ \ “ “k 8k 

k 

where it Is understood that in the present study repeated Indices imply 
summation, unless otherwise stated. Using Eqs, (2,11) and (2.12) the dis- 
placement field can be rewritten as 
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w = Jl x^ (3p + Bjj) + ^ n 

k ®k 



(2.15) 


Next Eqs, (2.13) through (2.15) are substituted into the combination of 

Equations (2.2) through (2.6). The resulting flap equation is multiplied 

? — 2 — 
by while the lag equation is multiplied by t dx^ and the equa- 

tions of motion are integrated over the domain 0 < x < 1. In this process 
the integrals associated with the elastic properties of the blade are 
integrated between 0 < x^ < £., while those associated with the aerodynamic 
loads are integrated between A < < -B. Thus the equations in their 

final nondimenslonalized form can be written as 

— ** __ *__2 f-s 

“fI 8l + 2 “Fi «Fi "sFj^ *1 + “f1 “f1 H ' [' *lk ^^Ikm ‘‘.J^k 

1=1, 2, ... M (2.16) 
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1 - 1, 2, ... » (2.17) 

where there is no summation over 1 on the left hand side of Equations 
(2.16) and (2.17).'^ 


In deriving Equations (2.16) and (2.17), the orthogonality relations 
for rotating beams given in Appendix A and the boundary conditions were 
used. The various coefficients used in the equations are generalized 
masses and are defined in Appendix B. The terms A^^ and A^^ are general- 
ised aerodynamic coefficients given by the expressions 


A 


Fi 



t|, dx 
1 o 


(2.18) 







(2.19) 


Next the aerodynamic terms will be specified* The loading terms in 
the 2-direction and the y-direction can be obtained from Friedmann 
(Ref. 13). 

= ap^bRU^ (U^e - Up) (2.20) 

r ^Do ”1 
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imr 



imr 


is usually zero. 
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where the velocities and are given by 

:[x ■ - 


Up = Ow + fiRj^X + V cos '{' j 


r 

X + y sin + y cos i(/ 


Sv 1 
9x J 

rs-* 


( 2 . 22 ) 


(2.23) 


The last term in Eq. (2.23) represents the component of tangential 
velocity due to radial flow along an elastically deformed blade. This 
term has been sometimes neglected even though it is not a higher order 
term. The additional terms due to this effect will be underlined by 
for future reference. 

For the purposes of this study* the inflow ratio X in Eq. (2.22) was 
evaluated using an expression for constant inflow ratio in hover . given 
by (Ref, 13) 

, aa I 
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1.^-1 

aa 


(2.24) 

This inflow relation is equivalent to taking the induced velocity 
at 3/4 blade radius as representative of a constant induced velocity over 
the whole disk. It is clear that for forward flight one should use the 
expression 


X “ y tan + C^/2 




A (2.25) 

However, use of this expression would have required the use of a trim pro- 
cedure by which 0 is changed as a function of y. This approach was used 
in Ref. 13, where it was concluded that the changes dictated by such a 
trim procedure have the tendency to mask the effects due to forward flight. 
Therefore, in order to illustrate clearly the effects of forward flight, 
it was decided that the use of Eq, (2,24) was more convenient. 
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By substituting Eqs. (2.13) through (2.15), (2.22) and (2.23) into 


Eqs. (2.20) arid (2.21), the expressions for and can be obtained. 
Carrying out the integrations in Eqs. (2.18) and (2.19) yields the general- 
ized aerodynamic loads and The complete expressions for Ap^ and 

1 24 

are given in Appendix C, The various flap coefficients F ,...,F and 
1 24 

lag coefficients L for the aerodynamic loads are defined in 

Appendix B. 

2.3 Linearization of the Equations 

The equations of motion will be linearized about a convenient equilib- 

f 

rium position which is chosen as the static equilibrium position in hover. 
The equations cannot be linearized about a forward flight condition having 
a particular value of p since there is no guarantee that such an equilibrium 
position exists as a point of stable equilibrium. 

The linearization is performed by writing 

o 


w^ = + Aw^ = 

^'^e = " [^m^m + 


(2.26) 


where the superscript denotes the static equilibrium position* 
From Eqs. (2.16), (2.17), (2.26) and Eqs. (C.3) and (C.4) in 
Appendix C 
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(l)= 


L? — + A9 - 
i a 1 i 


1 = 1 . 2 , 


.M 


(2.28) 


No summation over 1 on the left hand side of Eqs. (2.27) and (2.28). 


For every case stability of the equilibrium condition in hover is 
automatically checked. 


11 



Considering the case of one elastic mode for each degree of freedom 
and dropping the subscripts on the flap and lag cpefflclents, we have 




Eqs. (2.29) and (2.30) can be solved for the static 



II 
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equilibrium position 


where 
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\i 4.x - 
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(2.32) 


i("i" - Pd® + ^ (f)' ^ 

Using the order of magnitude analysis in Appendix C, only generalized 
aerodynamic terms up to and including second order terms are retained in 
the flap equation. Thus Eqs. (2.16), (2.31) and Eq . (C.3) in Appendix C 
become 
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From the order of magnitude analysis in Appendix C, only generalized 
aerodynamic terms up to and Including third order terms are kept in the lag 
equation. Therefore, from Eqs. (2.17), (2,31) and Eq. (C.4) in Appendix C, 


"l 1 % * P“Ll“il1sLl ^ "l (-2=111 + 2»V^ J - 26 b« 

■'■ 2^ (e) cosiKSp+ej,) f 9 + 1^3 2y si„|(, I ^ + 1,3* 2 I 
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+ ln,aL. - E® -t-X (Xr 


L uX cosUi6 
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2y cosip — 
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HP 


(1+cos 21^1) 
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(2.35) 


The third order terms in the lag equation are underlined by for 

future reference. Also, the radial flow terms are underscored with . 

Eqs. (2.34) and (2.35) are a system of coupled linear periodic equations 
under the influence of periodic forcing. The periodicity in the coefficients 
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and the forcing are due to the effects of forward flight. The equations are 
coupled through the aerodynamic, inertia and stiffness terms. 

Furthermore, for mathematical convenience the equations of motion 
have to be transformed into a system of first order equations. This is 
achieved by making the following substitutions, 
a|i - 

^ (2*36) 

“l ■ ^3 

Ahj - 74 

Thus, the equations of motion in their final form can be written as 

I = A W 2 + f (ip) (2.37) 

Where A is a 4 x 4 matrix and f is a column matrix. The elements 
of both are defined in Appendix C. 

The equations of motion (2.37) will have a different form for the 
normal and reversed flow regions. The representation of the reversed 
flow together with its effect on the form of Equation (2,37) is described 
in Appendix D, 
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SECTION 3 


MATHEMATICAL METHODS FOR DETERMINING 
THE STABILITY OF 
LINEAR PERIODIC SYSTEMS 

3.1 Introduction 

This section will consist of a brief summary of the methods available 
for obtaining the solutions of homogeneous systems of linear differential 
equations with periodic coefficients. In general such a system is governed 
by the following set of equations 

A 

y = A(iJ^) y (3.1) 

Where y is the state variable column matrix and A(i/^) is an n x n 
periodic matrix whose elements have a common period denoted by T» thus 

A(ij>+T) ^ kill)) (3,2) 

The problem of determining the stability of such a system has been 
considered in the literature associated with various fields: 

(a) mathematics (Refs* 14--16) 

(b) linear control system theory (Refs. 18-20) 

(c) dynamic stability or parametric excitation problems (Refs. 

21-25) 

(d) structural dynamics problems related to helicopter rotor blade 
dynamics (Refs. 1, 5, 6, 10^ 26-29) 

Unfortunately, there seems to be a considerable lack of communication 
between researchers working in these areas, which led to a considerable 
amount of overlap and duplication of efforts. Therefore it is believed 
that the various references mentioned above could be useful to other 
researchers dealing with similar problems. It should be mentioned that 
the list given above is far from complete and a complete review of the 
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bibliography dealing with equations with periodic coefficients is beyond 
the scope of the present study. From these references it is evident that 
in the past the following methods have been used: 

(a) Hill's method of infinite determinants (Refs. 21, 26 and 11) 

(b) Multivariable Floquet-Liapunov theory (Refs. 5, 6, 10, 18-20, 
23-25, 28) 


(c) Perturbation Methods (Refs, 8, 13, 21, 22) 

The first method has proved to be cumbersome and numerically ineffi- 
cient. The second method is the most general one, its main drawback has 
been the computational effort required for evaluating the transition ma- 
trix. This can be overcome by using Hsu’s method. The third method is 
limited to cases when the parametric excitation can be expressed in terms 
of some small parameter which tends to limit the generality of the method. 
3.2 Multivariable Floquet-Liapunov Theory 
3.2.1 The Transition Matrix and Its Properties 

The solutions of a system of differential equations with periodic 
coefficients are closely associated with the concept of the state transi- 
tion matrix, therefore it is important to start by defining this matrix. 

The transition matrix, or the fundamental matrix, is a matrix whose 
columns contain the linearly independent solutions to Equations (3.1). Thus 
the transition matrix in general is written as 




(3.3) 


where each of the columns in equations (3.3) satisfies Equation (3.1) 
with the initial conditions given by 


18 







j = k 

j k 


to yield • Thus clearly 


2<^o>V 



(3.4) 


(3.5) 


Two important properties of the transition matrix are given below 


(Ref. 18) 


det ) = exp 

o 



trace A (a) da 


(3.6) 


£(4».'I»q) = £(i|»,iJ^j^) i|>q) (3.7) 

3.2i2 The Floquet-Liapunov Theorem and Its Consequences 

Proof of the Floquet-Liapunov theorem is given in References 15 and 
18, for conciseness only the theorem Itself will be stated. 

Floquet-Liapunov Theorem . If A(i|M-T) = A(i|») then the transition 

»v 

matrix associated with Equation (3.1) can be written as 

^(<|J»i|^_) = P"^(i|))e?^'^“'*^o^P(i^ ) (3,8) 


where 

P(IHT) = P(i|;) (3.9) 

and R is a constant matrix. 

•V 

It can be shown that in general (Refs. 15 and 18) 

RT 

= C = J(T,0) (3.10) 

where C is also a constant matrix. Furthermore P(i|j) is given by 

P“^(<P) = $(i|;,0)e"5'^ (3.11) 

The most important consequences of this theorem are 
(a) Knowledge of the transition matrix over the period 
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0 < ij; < T determines the solution to Equation (3.1) everywhere 
because from Equations (3.7) and (3.10) 

£(<|H-T,0) = (3.12) 

and in general 

( ' RT\^ 

e~ ) (3,13) 

where s is any integer and 0 < ^ < T, 

Therefore using Equation (3.11) P(i|;) is determined everywhere, 

(b) Knowledge of the transition matrix at the end of the period determines 
the stability of the system. From Equation (3.10) the transition 
matrix at i/; = T clearly determines both C and R. Here two cases need 
to be considered. 

Case (1) The matrix C had n-independent eigenvectors associated with 
n-distinct eigenvalues. For this case this also means that R has 
n-independent eigenvectors, therefore from elementary linear algebra a 
similarity transformation can be found such that 

9 "^ 5 2 “ i (3.14) 

where the columns of Q are the n-linearly independent eigenvectors of R and 
A is a diagonal matrix whose elements are the eigenvalues of R. Combining 
equations (3,10) and (3.14) and using the definition of the matrix exponen- 
tial (Ref. 15) one has 



or 

- A » q“^ C Q * q"^ $ (T,0)Q (3.15) 

where A is a diagonal matrix containing the eigenvalues of the transition 
matrix at the end of one period. The eigenvalues of $ (T,0) or the charac- 
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teristic multipliers are related to the eigenvalues of R, denoted char- 


acteristic exponents ^ through the relation 


V 

e = k = 1, 2, ... n (3.16) 

Clearly and are both complex quantities in general, thus 

(3.U) 


from which 





(3.18) 

(3.19) 


the quantity can be determined according to the Floquet-Liapunov 

theory only within an integer multiple of the nondlmensional period. 

The stability criteria for the system is related to the eigenvalues 

of R or the real part of the characteristic exponents C . The solutions 

k 

of the Equation (3.1) approach zero as if 


I 4r "kl I ^ 1 


< 0 k . I, 2, 


n (3.20) 


Case (2) The eigenvalues are not distinct anj^ore. In this case instead 
of finding a similarity transformation which transforms R into diagonal 
form, a similarity transformation can be used which transforms R into the 
Jordan canonical form, this case is treated in Reference 15. 

3.3 Methods for Calculating the Transition Matrix 
3.3.1 General 

From the description given in the previous sections it is clear that 
in order to obtain the solutions for a periodic system using the Floquet- 
Liapunov theorem one has to evaluate the transition matrix over one period. 
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This can be done using the methods described In the following sections, 

3.3.2 Direct Numerical Integration 

This straightforward method, which has been used in Refs. 6 and 28 

is based upon equations (3.3) and (3,4), The various columns of the 
transition matrix are evaluated by direct numerical integration of the 
equations of motion, Equations (3.1), using a fourth order Runge-Kutta 
method, with the initial conditions specified by Equations (3.4), 

3.3.3 Hsu's Method for Evaluating the Transition Matrix 

In a recent series of papers, Hsu (Refs, 23-25) has developed various 
methods for approximating the transition matrix during one period; the 
most efficient one consists of approximating the periodic matrix A(4)) by 
a series of step functions, this method can be considered to be the gen- 
eralization of the "rectangular ripple" (Ref. 1) method to multi- 
dimensional systems. The method consists of evaluating the state transi- 
tion matrix by dividing a period into a number of equal parts and con- 
sidering the equations over each interval to be a set of constant coef- 
ficient equations. Due to the fact that the method is new and numer- 
ically efficient a description of its essential aspects is given below. 
Each period T is divided into K intervals denoted by 
k = 0, 1, 2....K, with 0 = < 

The k~ interval is denoted by and its size by 

“ '^k-l* — interval the periodic coefficient matrix 

A(iii) is replaced by a constant matrix which is defined by 

^kC \ (3.21) 

1 Pk 

^k = r X 

^ \ ^k-1 " 
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Thus the actual system, Eq. (3.1), is approximated by an approximate system 
x(tjj,k) = C(\jj;k) x(i|;;k) (3.22) 

where 

x(0;k) “ C 

and 

oo k 

C(t|^;k)=>X) Z Cj^[U(4»-sT-i|;j^_^)-U(i{/-sT-if/ )] (3.23) 

^ 111= k“l 

The elementary theory of differential equations with constant coef- 
ficients enables one to write the fundamental matrix (or the transition 
matrix) of the system with $.(0,k) = I as 

= exp [(V^-’^K_i)CKjexp(^_^Aj,_j^)....exp(C^^^ (3.24) 

and the growth matrix, or approximate transition matrix at = T is given by 

K 

H(K) = exp(A^C^) e*P j^) ... exp(Aj^Cj^)_ = IT exp(A^C^) (3.25) 

~ i=l ■" 

With regard to the product sign, it is understood that the order of 

positioning of the factors is material and that the k— factor is to be 

til 

placed in front of the (k-1)— factor. 

It is shown in Reference 25 that when K -»■ «> 


x(lf^,K) y(M>) 

•Nr ~ 

2a (^.K) ^ $(!{», o) 

and 


H(K). 


H(K) = ^ 

The basic numerical problem is therefore the efficient computation of 

» 

Using the definition of the matrix exponential 


(A C )" (A C ) 

exp(A^C^) = I + A^C^ + — + ... 


n 

> -h « • < 


(3.26) 


XI CO 
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For small time intervals o and the series in Equation (3.26) 

converges rapidly, and the value of the matrix exponential ..can he accu-_ 


rately approximated by a finite number of terms. Thus 

^ (AfC “ (A C 

exp(A,Cl^I+E -T^ =exp(AC)-2^ - 

~ ^=1 i- j=j+l J! 


Therefore the approximate value of the growth matrix can be written as 
K ( -J (A.C,)^ 

H ,(K) = n I + £ v7^ . (3.28) 

i=l j=l 

General error bounds for these approximations are obtained in Ref. 25, 
furthermore it can be shown that 

$(T,o) - H^(K) = 0 (A^) (3.29) 

for J > 2. 
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SECTION 4 


RESULTS AND DISCUSSION 

4.1 Numerical Quantities Used in the Calculations 

In computing the numerical results the following assumptions were 

made: 

Mass and stiffness distribution was assumed to be constant along 
the span of the blade. 

Two different kinds of mode shapes were used: 

(a) For the cases when essentially only a study of the trends or 
properties of the numerical methods (i^e. numerical experimenta- 

. tion) was conducted, the mode shapes in flap and lag were both 

approximated by the first nonrotating assumed mode shape (Ref. 13) 

n, = Yi “ - (l-x )^] (4.1) 

11 3 o o 

This equation satisfies all the boundary conditions of the 
problem. When an assumed mode shape is used the elastic coupling 
effect as represented by Equations (E.5) through (E.8) was neglected. 

(b) For the majority of the cases, where a parametric investigation 
of typical hingeless rotor blade configurations was conducted, an 
exact rotating mode shape in flap and lag was employed. The exact 
rotating mode shapes were generated by using Galerkin^s method 
based upon five nonrotating cantilever mode shapes for each flap or 
lag degree of freedom. For these cases the elastic coupling effect 
was included. The method by which the rotating modes were obtained 
is essentially identical to that described in Reference 31 and the 
various numerical constants required were taken from Reference 32, 
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The flap and lag coefficients F^, and the various generalized 
masses defined in Appendix B were calculated using seven-point Gaussian 
integration. The numerical values for these coefficients were calculated 
for i/R » 1.0, A = 0.0, B = 1.0 and e^^ = 0. 

As mentioned before, the inflow was evaluated using Eq. (2.24). It 
has been shown in Reference 13 that the correct procedure would involve a 
trim requirement for which 9 must be changed as a function of p at a 
fixed value of C^. This added complication was not considered to be worth- 
while in the present trend type study. 

Finally, in all the computations the following values were used; 

= .01; a = 2Tr; a = 0.05 

Various other pertinent quantities are specified on the plots, 

4,2 Brief Description of the Computer Program 

This section describes briefly the computational implementation of 
the mathematical techniques described in Section 3. As pointed out in 
Section 3, the method of solution is based upon the Floquet-Liapunov 
theorem, but the transition matrix at the end of a period is evaluated 
using two different approaches: (a) Direct numerical integration, 

(b) Hsu*s method of approximating a periodic system by a series of 
step functions. Consequently, the two computer programs by which the 
results were obtained are identical except for the routine which evaluates 
the transition matrix. 

The computational steps are outlined below: 

1. The generalized masses, and generalized aerodynamic coefficients 
i 1 

F , L are evaluated using seven point Guassian integration. In 
the calculation of the aerodynamic coefficients the normal and 
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reversed flow region are accounted for as described in 
Appendix D. 

The transition matrix at the end of one period is computed using 
two separate methods 

(a) Direct Numerical Integration. For this case the elements 

of the A-matrix are generated and the resulting equations are 
integrated using the fourth order Runge— Kutta method available 
in the IBM~Scientif ic Subroutine Package. The numerical inte- 
gration is performed four times yielding the four columns of 

the transition matrix. The time steps used were the same as 

/ 

those used for Hsu^s method. 

(b) Hsu’s Method. The transition matrix is evaluated using Equa- 
tion (3.28). After some numerical experimentation it was found 
that A = T/50 gives good accuracy. The periodic part in the 
elements of the A-matrix (Appendix C) were evaluated using 
Equation (3.21) which yields 

k— 1 




A 

f 


cos d i/j = cos m 


kavg MiJ) 


■ fM\ 

sin \^\ 


k-1 


r 


sin 2 ^ d = sin ^ sin Aij; 


k-1 


_A r 

Wi\_^ 


cos 2 d = cos ^ sin 


where 




*k + I'k-i 


kavg 


: At - 


(4.3) 


(4.4) 

(4.5) 
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It is apparent from these equations that for small intervals 
the periodic terms could be replaced by their value at the 
mid-point of the Interval, 

3. The eigenvalues of the transition matrix at the end of a period 

are evaluated using a Jacobi-type eigenvalue routine. The charac- 
teristic exponents are then evaluated from Equations (3.18) and 
(3.19). Of the four characteristic exponents obtained, two are 
associated with the flap degree of freedom, and two are related to the 
lead-lag degree of freedom. The critical advance ratio is the 
advance ratio for which one of becomes positive. The degrees 

of freedom are identified by tracing back their history from y = 0, 

4,3 Results 

The results obtained from the present study are usually given in 
the form of plots representing the variation of the real part of the 
characteristic exponent with the advance ratio y. 

Before proceeding with the description of the results it is important 
to note that Figures 3 through 8 are based upon the assumed mode shapes 
Eq, (4.1), and in these cases the elastic coupling effect is neglected. 
Figures 10 through 17 have been obtained using the exact rotating mode 
shapes and in these cases the elastic coupling effect is included. 

Figure 3 shows results for a typical case using Runge-Kutta numer- 
ical integration and Hsu’s method with different orders of approximation 
used in calculating the state transition matrix. This plot shows that 
Runge-Kutta and Hsu's method using a fourth order approximation to the 
solution of the constant coefficient equations in each interval give 
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almost Identical results « As shown, the lag degree of freedom becomes 

unstable for y = .754 and the frequency of oscillation is w, = 1,286. 

k 

This result clearly indicates that by neglecting the lag degree of 

freedom one could obtain completely incorrect stability boundaries. 

Figure 4 is a plot of the C.P.U. time required on an IBM 360/91 

computer to calculate the Floquet transition matrices for the four cases 

depicted in Fig. 3. As shown on the figure, Hsu's method, even using a 

fourth order approximation, is computationally much more efficient than 

the fourth order Runge~Kutta numerical integration scheme. 

Figure 5 is a comparison of results between fourth order Runge- 

Kutta and Hsu's method for a case that has a relatively high value of u . 

c 

Again the results compare favorably even for high values of y. 

In the derivation of the equations of motion, some new terms due to 

the effect of radial flow were included,* 

Figures 6 and 7 show the effects of the radial flow terms on the 

stability of the system. The effect of neglecting the radial flow terms 

on the real part of the characteristic exponent associated with the flap 

degree of freedom is illustrated in Fig, 6. Note that the radial flow 

terms have a stabilizing effect on the flapping motion; with the radial 

flow terms suppressed, the flap degree of freedom becomes unstable 

at y = 1.33. The effect of the radial flow terms on the lead-lag degree 

of freedom is illustrated by Fig, 7; as shown, without the radial flow 

terms the instability is completely eliminated and the system becomes 

unstable in flap. When the radial flow terms are included, the lag degree 

of freedom is the critical one and the system becomes unstable at y = .754, 

*The radial flow terms in this discussion are those due to the under- 
lined term in Eq. (2.23). 
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This matter was pursued further to determine which radial flow term was 
causing the instability in the lead-lag degree of freedom. The actual 
destabilizing term in the equations of motion was found to be an aerody- 


namic coupling term. This term couples the flap motion with the lag 


motion in the flap equation. Its form is 

9w 9v 
8x 9x 


2 2 3w av 

]1 cos ^ ^ 


(4.6) 


o o 

this term is due to the term in Equation (2.21). The term shown 

above is the complete nonlinear one, clearly the one retained in the 
equations of motion is the coupling term obtained from linearizing this 
expression. 

It has been pointed out by Dr. K.H. Hohenemser, at a recent meeting,* 
that this term could be transformed into a constant coefficient term, if the 
equations of motion are rewritten in a nonrotating coordinate system such as 
multiblade coordinates. 

3 

As outlined in Appendix C, the aerodynamic terms associated 

with the lag equation were retained in the derivation. Fig. 8 shows the 
effect of these third order terms on the real part of the characteristic 
exponent associated with the lag degree of freedom plotted as a function 
of the advance ratio y. Clearly the third order terms associated with 
damping are the relatively important ones. The third order terms asso- 
ciated with damping in the lag equations are underlined by in 

element A^« of A, Appendix G, 

Figure 9 is given for convenience. It yields the relation between 
the first rotating and nonrotating nondimens ional flap and lag frequen- 
cies, respectively. 

_ 


AHS /NASA-Ames Specialists Meeting on Rotorcraft Dynamics, Moffett Field, 
Calif,, Feb. 13-15, 1974. 
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As mentioned, the results presented in Figs. 3 through 8 were 
obtained by using the assumed mode shape and neglecting the elastic 

coupling effect. In order to assess the effect of this assumption the typical 
case has also been recomputed with the exact rotating mode shape and the 
elastic oupling effect, the results are shown in Figure 10* From Fig. 10 it 
is clear that use of the exact rotating mode in flap and lag reduced the 
value of y to y = ,653, when the elastic coupling is also included y is 

C C IT O ^ 

further reduced to y^ = *583. Thus, for this case, y^ seems to be more 

sensitive to the type of mode shape used than to the inclusion of the elastic 

coupling effect. It is also interesting to note, that for this case the 

elastic coupling effect is destabilizing, while for hover its effect on 6 is 

c 

quite stabilizing. (Ref. 30) 

It should be noted that all the results starting from Fig. 10 are 
based upon the exact rotating mode shapes and the elastic coupling effect 
is always included. 

The importance of the reversed flow region is illustrated by Fig. 11. 

As shown, with the reversed flow region the instability occurs at a higher 
value of y than without the reversed flow region. Similar trends were 
observed in previous studies when only the flapping motion was considered 
(Ref. 5), indicating that by neglecting the reversed flow region one could 
expect conservative results from a stability point of view. It also 
indicates that for this particular case the effects of reversed flow are 
negligible for y < 0.6 and start being important for y > 0.8. 

It is also clear from Fig. 11 that the difference in due to 
inclusion of the reversed flow region, for y = 1.0, is approximately 7%. 

Thus, the approximate representation of the reversed flow region, described 
in Appendix D, cannot significantly affect the accuracy of the results 
presented in this study. 
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Previous studies (Ref, 12) dealing with the effect of viscous type 
of structural damping on the stability boundaries- for hover indicated 

this parameter has an important stabilizing effect. The effect of 
this parameter for forward flight is shown in Fig, 12. The stabilizing 
effect of structural damping In the lag degree of freedom is evident 
from Fig. 12 where the value of y^, at which the characteristic exponent 
associated with the lag degree of freedom becomes positive, is plotted 
as a function of the viscous type of structural damping in lag. It is 
interesting to note that this plot indicates that the greatest stabilizing 
effect due to damping is obtained in the range 0 < < .02 (2% of 

critical damping), after which the gain in stability tends to level off. 
Similar trends were also obtained from stability studies in hover. 

Again in order to Illustrate the sensitivity of the results to the 
mode shape used and the elastic coupling effect three sets of results 
are presented. As seen the results are mode shape dependent and inclusion 
of the elastic coupling effect tends to reduce the value of y^, at which 
the instability occurs. 

The effect of collective pitch setting on blade stability is illus- 
trated by Fig. 13 for the typical case (0^^^ = 1.175 and = 1,28303, As 
shown for a low value of collective pitch setting 9 = 0,05 the lag degree 
of freedom becomes unstable at y^ == 1.435 and “ 1,279. For a high 
value of collective pitch setting 0 = .30 the lag degree of freedom becomes 
unstable at y^ - .535 with a frequency of = 1.332, Similar results 
have been obtained for a variety of other cases, the results are not 
given here for the sake of conciseness. These results consistently 
indicate that there is a degradation in stability in forward flight with 
the increase in collective pitch setting. These results also Indicate 
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that the use of a trim procedure, such as developed in Reference 13, 
should be required in realistic stability calculations in forward flight. 
Furthermore, it seems that the assumption of nonlifting rotors used in 
some forward flight studies (Ref. 7) can be nonconservative. 

Figure 14 illustrates the effect of the Lock number, based on normal 
flow, on the stability of the blade in forward flight. As shown, an 
Increase in the value of y tends to destabilize the blade in the lag 
degree of freedom. For y - 5^0 the lag degree of freedom becomes unstable 
at = 1.425 and = 1.283, for y = 10. the lag degree of freedom 
becomes unstable at = .583 and = 1.295, and for y = 15. the lag 
degree of freedom becomes unstable at y = 0,568 and o), = 1^294. 

Hingeless helicopter blades are usually designed with a first rotating 
flap frequency of 1.1 ^ l.#2. Thus the value of = 1.175 used in 

the numerical calculations can be considered to a representative value 
for this parameter. On the other hand the rotating lag frequency for 
hingeless blades has usually a considerable variation, therefore it is 
reasonable to investigate the effect of changing this parameter over a 
wide range such as 0,2 ^ _< 2.5, the results are presented in Fig, 15. 

For the sake of completeness two sets of results are presented, those 
with the assumed mode shapes and those with the exact rotating mode 
shapes and the elastic coupling effect. 

Two types of instabilities are shown in Fig. 15: (a) Those for 

which the imaginary part of the characteristic frequency is equal to 0, 

1/2 or 1 indicating that the instability is directly due to the periodic 
coefficients in the system, for these points the appropriate number 
appears inside the circles or squares indicating the flutter points and 
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(b) Those for which o)j^ has a continuous value, usually close to the 
rotating lag frequency, marked by letters A, B, C etc. 

As can be seen from Fig. 15 there are two regions where degradation 

of stability can be expected with forward flight. One region is between 
_ 

0.9 ^ ^1.6, while a second and much smaller region is between 

0,45 < < 0.55. The detrimental effect of forward flight for this 

region has been predicted also from purely mathematical considerations 
(Ref. 9). 

The results presented in Fig. 15 seem to indicate that at high 
advance ratios a soft inplane hingeless blade with 0.7 < < ,8 will be 

the most stable; while for a stiff inplane hingeless blade a lag frequency 
of > 1,6 would be the best. 

Finally, it should be noted that, although the points in Fig. 15 
have been connected by lines, a possibility exists that additional flutter 
points which could be evaluated between the calculated points may not be 
exactly on the lines as drawn, this is due to the periodicity of the coef- 
ficients in the equations of motion. 

Previous studies have indicated that droop and preconing can 
significantly affect the value of 0^ at which the linear system becomes 
unstable in hover. It was reasonable therefore to investigate how these 
parameters effect the stability of the blade in forward flight. 

The results for the typical case are shown in Fig. 16, for this 
particular combination of flap and lag frequencies is relatively 
insensitive to variations in 8p, while negative angles of droop seem 
to be stabilizing and positive angles of droop seem to be quite 
destabilizing. 
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The effect of the same parameters on a case with a rotating lag 
frequency close to ~ 1»0 is shown in Fig. 17. In this case both 
droop and preconing are strongly destabilizing. Thus it is clear that 
the effect of these parameters is strongly related to the combination of 
rotating flap and lag frequencies for the particular blade configuration. 



SECTION 5 


CONCLUSIONS 


The most Important conclusions which can be drawn from this study 

( 

are summarized below: 

(1) Hsu^s method for calculating the Floquet transition matrix at the 
end of a period is a powerful computational scheme which represents 
a significant improvement over previous methods. Compared to direct 
numerical integration, Hsu^s method was approximately 2*5 times 
faster for a system of four first order equations. For higher order 
systems, it probably could be more efficient by one order of 
magnitude » 

(2) Flapping instability and response studies at high advance ratios can 
be inaccurate and misleading due to the neglect of the lead-lag 
degree of freedom. The effect of the periodic coefficients on the 
coupled flap-lag system shows that, in general, instability can 
occur at lower values of advance ratios than when the flap degree of 
freedom is considered by itself. 

(3) In addition to the instabilities associated'with the periodic 
coefficients (l.e. with frequencies of 0, 1/2 or 1), the coupled 
flap-lag system has the tendency to become unstable due to an aero- 
dynamic coupling effect associated with the radial flow terms. This 
instability, which has a frequency close to the rotating lag frequency 
of the system, usually occurs at values of much lower than those 
for which the flapping degree of freedom becomes unstable. 

(4) The effect of the reversed flow region is negligible for p < hut 

it becomes quite important for y > 0*8. The approximation of the 
reversed flow region, outlined in Appendix D, seems to be acceptable. 
The approximation is better for high advance ratios. 
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(5) Viscous type o£ structural damping in the lead-lag degree of freedom 
has a strong stabilizing effect on the instability due to the aerody- 
namic coupling effect associated vlth the radial flow terms. 

(6) The value of the collective pitch setting has a considerable effect 
on the value of the advance ratio at which instabilities due to the 
periodic coefficients or the radial flow aerodynamic coupling terms 
occur. An Increase in collective pitch is destabilizing, therefore 
high advance ratio studies which do not include this effect 
(nonlifting rotors) may yield unconservative results. 

(7) The blade Lock number has a destabilizing effect on flap-lag 
instability due to aerodynamic coupling effect associated with 
radial flow terms. 

(8) The results obtained seem to indicate that certain values of rotating 
lag frequencies can provide improved aeroelastic performance for 
hingeless blades at high advance ratios. 

(9) Droop and preconing can significantly affect blade stability in 
forward flight. The effect of these parameters is dependent upon 
the particular flap and lag frequencies of the blade configuration 
being considered. 
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Figure 1. Displacement Field of a Torsionally Rigid Cantilevered Blade 
with Droop and Preconing. 




fk CHARACTERISTIC EXPONENT FOR LAG 



Figure 3. Comparison of Runge-Kutta and Hsu's Methods for a Typical Case. 


CPU TIME TO COMPUTE TRANSITION MATRIX (SEC) 



Figure 4. Variation of CPU Time to Compute Transition Matrix With Advance Ratio. 
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Figure 5. Comparison of Runge-Kutta and Hsu's Methods for a Case With High 
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CHARACTERISTIC EXPONENT 






CHARACTERISTIC EXPONENT 


1. ALL 0 (eg) TERMS IN LAG EQ. NEGLECTED 

2. ALL TERMS INCLUDED 

3. ONLY 0 (eg) TERMS ASSOCIATED WITH 



Figure 8. Effect of Third Order Terms in the Lag Equation on Characteristic Exponent for Lag. 
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CHARACTE-RISTIC EXPOIMEWT 



Figure 10. Effect of Exact Mode Shape and Elastic Coupling on Characteristic Exponent for Lag. 
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- CHARACTERISTIC EXPONENT IN LAG 



Figure 11. Effect of Reversed Flow on the Characteristic Exponent for Lag. 


c3p,= 1.175 
ool^= 1.28303 



Figure 12. Effect of Viscous Structural Damping on 
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- CHARACTERISTIC EXPONENT IN LAG 



Figure 13. Effect of Collective Pitch Setting on Typical Case. 
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Figure 15. Effect of Variation in Rotating Lag Frequency on n . 

c 
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Figure 16. Effect of Droop and Preeoning on for Typical Case. 
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Figure 17. Effect of Droop and Preconing on 
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APPENDIX A 


ORTHOGONALITY CONDITIONS 

The orthogonality condition for a rotating beam vibrating out of its 
plane of rotation is obtained from Friedmann (Ref, 13). with 0=0 
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The orthogonality condition for a rotating beam vibrating in its plane 
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APPENDIX B 


preceding 


COEFFICIENTS ASSOCIATED WITH GENERALIZED MASSES 
DAMPING AND AERODYNAMIC LOADING 
The non-dimens Iona li zed general masses are 
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The viscous type of structural damping ratio is defined by 
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Listed below are the coefficients associated with aerodynamic loads. 

These integrals are evaluated between the lower limit A and the upper limit 
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B where A and B are the tip loss factors. The calculation of the aerodynamic 

loads in the mixed or reversed flow regions is described in Appendix D. 
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appendix c 


GENERALIZED AERODYNAMIC LOADS AND ELEMENTS 
OF THE COEFFICIENT MATRIX AND FORCING VECTOR 
C.l Generalized Aerodynamic Loads 

In evaliiating the generalized aerodynamic loading terms, one encounters 
a considerable number of terms which are small and therefore negligible. In 
order to neglect the appropriate terms, a rational ordering scheme is used 
which enables one to neglect terms in a systematic manner. In this scheme 
all the Important parameters of the problem are assigned orders of magnitude 
in terms of a typical displacement quantity thus: 
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An order of magnitude analysis of the equations indicates that in 

general terms up to and includiiig ) must be included in the linearized 

flap equation, while for the lag equation some terms have to be 

retained. These higher order terms in the aerodynamic part of the flap and 

lag equations will be marked by a y^and neglected. 

The aerodynamic loading terms L and L are given by the expressions: 
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Using these equations and the definition of the flap and lag coefficients 
given in Appendix B, the generalized aerodynamic forces are 


^1 


^ -iij 


fJ 0 + 2y 6 sini/j 


,3 y 


- F‘ y cosi|j(Bp+Bjj) + ^ 0(l-cos 2ip) “ sin 2tKBp+Bp) 


- X - F^ Xy sin?{) + - F^j^ y cosi|> - F^^^ ^ sin 2\jj 

+ F- F® i- F^ 

L ik R ik 


Sk 


^ P 2T. 22 

y — sini^gj^ + - F^^ 2y cosijj0+F^^ yX cos\li 


im 


- F-f 6 sin 2^ + F?! ^ ce^+ej (1+cos h^ + [- fJ 


xm 


,11 ^ 

R 

2 


_im 1 P 

,11 ' Z 


“ 2|e 

xm R 


,11 


- ^im ^ R 2 slntl^e + y - cos^(8p+Bj,) + F^ X 


il h 
rJ 


* h f U cos* 


,24 


Z * 


15 IZ\^ * * 


^ ^ll^m ^ i ®k "m + t ^k \ 


(C.4) 


69 



“ ■ “b - A ■5' 


X V ^1 ■'*<> ■ 2 ('■i ■" ^1 >‘< 6 p + eD >'^'>®*0 


Sc 


C 2 

. _ 1 - Do .* ^2^2 ^2li . „ _. ,_ .-? 0 M 

+ 2y — sinyj - a ~ ^i t~ ^^p'''*^D^~ '*' f~ (1-cos 2\p) 

2 

- L^ A 2ii(Bp+6jj)cos4' + L^ Xu© slnif) + L^ (gp+3p)6 sin 2iJj 


i a 


2 

U cosi|;0 + ^ sin 2i|;0 - y^(0p+6_) (l+cos 2^) 


- 2Xu cosi|j 


ik ' “ik 2 ■ “Ik ^ 

7 a ^ . ,8 I . . ,8„5, 




^ik R ® •*• ^Ik R ^ - ^ik R 


- T ® 2 X ^ 
^ik R 


gk + 


Xu cosij^8 - Ljm 2~ (^^13)6(1+^08 2\p) 


^20 ,2 So , _, _23 „ , So"^ 

XI — T" 9in 2tP - L. 2u cosi|j 

im a Im • ^ a 


h + 
m 


- u co»K ye„) f e 


- lJ^ 2u slni() -| — - 2 f — - lP A 4 9 

tm ^ R a im R a im R 




m 


^ j o ^ 

■ ^iU 2^ 2# gj^gjj - ^ (1+cos 2ij;)eg^h^ - ^Ik£ 2 ^ ^ kh 


.16 


22 a 


■ Skm ^ P ® 8t ^n, " F V cosif)0 g,^ h. 


R " ^k • m "ikm R ^k m 


..18 /n * * 

ikil \R/ ®k 


- ("aV 

ikm I R/ 


* * 

0 81. 

®k m 


CC.5) 


C. 2 Elements of the A Matrix 

The elements of the A matrix, which defines the equations of motion when 
written as first order differential equations, are given below 
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C .3 Elements of the f Vector 

The elements of the f vector, which defines the forcing terms when the 
equations of motion are written as first order differential equations, are 
given below 
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APPENDIX D 


REGION OF REVERSED FLOW 


The circular region of reversed flow, which exists over the retreating 
blade, is quite well known. In past treatments of reversed flow it has been 
customary (Ref. 3) to define three separate regions: (a) norrmal flow, (b) 

(b) reversed flow, (c) mixed flow, and evaluate the appropriate aerodynamic 
expressions for each region. When this model is used together with a modal 
representation of the blade, the evaluation of the generalized aerodynamic 
expressions F^, becomes quite cumbersome, and a more convenient procedure 
had to be devised. 

The approximate reversed flow model developed for the present study 
consists of replacing the circular region by an approximate region which has 
two straight boundaries and a circular one as shown in Fig. 2, The approxima- 
tion is based on the assumption that the area contained in the circular region 
must be equal to the area contained in the approximate region. Two separate 
cases must be considered: (1) p < 1, (2) y > 1. 

Case (1). For y < 1, the diameter of the circular sector is taken as y. 

Simple geometric considerations show that the angle a is always a 
constant given by 

a = Tj/2 (D.l) 

Case (2) . For y > 1, simple geometric considerations show that 

a = 7T - 2 sin“^|^j + y^ sin"^ (D.2) 

Thus for li < 1 the generalized aerodynamic loads in the mixed flow region 
are calculated from ' 
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These equations are based on the assumption that the lift curve slope in 
the reversed flow region Is equal to the negative value of the lift curve 
slope in normal flow. 
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APPENDIX E 


ELASTIC COUPLING 


The angle of collective pitch 6 is the source of the so called elastic 
coupling effect in cantilevered rotor blades. This effect couples the bending 
perpendicular to the hub plane with the bending parallel to the hub plane. 

For convenient representation of this effect* the elastic coupling terms in 
Eq. (2.1) must be rewritten 


(EI)y 

cos^e + (El) sin^e = (El) + KEI) -(El) 1 sin^0 
z y L z 'yj 

(E.l) 

(EI)y 

sln^e + (El) cos^e = (El) - f(EI) -(El) 1 sln^0 

z z L z yj 

(E.2) 

Def inlng 



"cl“ 

[(EI)^-(EI)J sin^e 

(E.3) 

®c2 = 

j^(EI)^-(EI)J sine cose 

(E.4) 


When applying Galerkin*s method on the equations of motion, it is 
convenient to define the following additional expressions 
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